Spontaneous Crystallization And Filamentation Of Solitons In Dipolar Condensates 
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Inter-site interactions play a crucial role in polar gases in optical lattices even in the absence 
of hopping. We show that due to these long-range interactions a destabilized stack of quasi-one 
dimensional Bose-Einstein condensates develops a correlated modulational instability in the non- 
overlapping sites. Interestingly, this density pattern may evolve spontaneously into soliton filaments 
or into a checkerboard soliton crystal that can be so created for the first time in ultra-cold gases. 
These self- assembled structures may be observed under realistic conditions within current experi- 
mental feasibilities. 
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I. INTRODUCTION 

Recent experiments are opening new avenues for 
the study of the fascinating physics of dipolar gases 
[1, 2]. These gases present a significant electric or 
magnetic dipole-dipole interactions, which being long- 
range and anisotropic differ significantly from the short- 
range isotropic interactions usually dominant in quantum 
gases. Ultra-cold polar gases in optical lattices are par- 
ticularly interesting. Contrary to the non-dipolar case, 
polar lattice gases are characterized by significant non- 
local inter-site interactions that result in a rich variety of 
novel physical phenomena [2, 3]. Remarkably, the inter- 
site interactions play a crucial role even in deep lattices 
where hopping is negligible. In particular, dipolar Bose- 
Einstein condensates (BECs) in non-overlapping lattice 
sites share common excitations modes. This collective 
character enhances roton-like features in the excitation 
spectrum [4] and modifies the BEC stability, as recently 
shown experimentally [5]. 

Quasi- ID geometries allow for the existence of BEC 
solitons and hence modulational instability in these sys- 
tems leads to the formation of ID patterns, so-called 
soliton trains [6]. On the contrary, dynamical instabil- 
ity in higher-dimensional BECs is typically followed by 
condensate collapse [7]. In consequence, solitons patterns 
in higher dimensions, as e.g. a 2D crystal of solitons, are 
fundamentally prevented in non-polar BECs. 

In this paper we show that the destabilization of 
a dipolar BEC confined in a stack of non-overlapping 
quasi- ID tubes may be followed by the spontaneous self- 
assembly of stable soliton filaments or a 2D checkerboard 
crystal of solitons, providing a route for the first realiza- 
tion of self-sustained 2D arrangements of BEC solitons. 
This dynamical self-assembly stems from the correlated 
character of the corresponding modulational instability. 
While for non-dipolar condensates the instability in each 
lattice site would develop independently, the non-local 
dipolar interactions couple the non-overlapping BECs to 



form a density pattern shared among all sites. As we 
show, correlated modulational instability may be observ- 
able in current Chromium [8] and Dysprosium [9] exper- 
iments. 

The dynamically formed soliton filaments resemble 
dipolar chains of classical dipoles [10], as well as chains 
predicted for polar molecules [11, 12]. However, com- 
pared to the latter, soliton filamentation is expected to 
occur for smaller dipole moments due to the many-body 
character of each soliton. Remarkably, inverting the sign 
of the dipolar interactions results in the development 
of an anti-correlated density pattern which may be fol- 
lowed by the spontaneous formation of a stable crystal 
of solitons. This 2D checkerboard crystal resembles the 
Wigner-like crystal predicted for polar molecules [13, 14]. 
However, contrary to the latter, it is dynamically formed 
and self-maintained by a non-trivial interplay between 
intra-tube attractive and inter-tube repulsive dipolar in- 
teractions. 



II. MODEL 

We study below a dipolar BEC confined in a stack 
of quasi-lD tubes formed by an optical lattice (Fig. 1). 
The lattice is assumed to be sufficiently deep to suppress 
inter-site hopping. In each of the N m lattice sites the 
^-confinement is approximated by a harmonic poten- 
tial with frequency uj ± , whereas for simplicity we assume 
no confinement along z direction. We consider atoms 
with a magnetic dipole moment /i (the results are equally 
valid for electric dipoles, as e.g. polar molecules) ori- 
ented along y direction by an external magnetic field. 
The dipoles interact with each other via the dipole-dipole 
potential Vd (r — r') = g d (l — 3 cos 2 6) / |r — r'| , where 
gd = /ioM 2 /47r, with //q being the vacuum permeability 
and the angle formed by the vector joining the two 
interacting particles and the dipole moment direction. 

We assume the chemical potential much smaller than 




FIG. 1. (Color online) Scheme of the stack of disjoint quasi- 
1D dipolar BECs. 



Jtlj ± (this assumption is self-consistently verified in our 
calculations). Hence, we can factorize the BEC wave 
function at each site j, \£j (r) = (j)j (x, y) ipj (z), with 
(j)j (x, y) the ground-state wave function of the xy har- 
monic oscillator. Treating the dipolar potential in the 
Fourier space [15] we arrive at a system of N m coupled 
ID Gross-Pitaevskii equations describing the BEC stack: 
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l ± = yJh/muj_ L is the xy oscillator length, A is the lattice 
spacing, and g = Aiiah 2 /ra. Note that for a fixed ratio 
A/Zjl the physics of the system is governed by the values 
of g and g d . 



III. LINEAR REGIME: BOGOLIUBOV MODES 

Starting from a homogeneous on-site linear density no 
we are interested in the dynamics that follows the desta- 
bilization of the condensate after an abrupt change of the 
scattering length a. A substantial insight into the first 
stages of the post-instability dynamics is provided by the 
analysis of the elementary excitations of the condensate. 
To this end we introduce a perturbation of the homo- 
geneous solution, ipj (z,t) = [y/riQ + \j ( z ^)] e~ %fIjt / h , 
with Xj(z,i) = Uje< z «- Ut ) + v * e -< z «- ut ) , where fij is 
the chemical potential in a site j, and q and uj are the 
z-momentum and the frequency of the elementary exci- 
tations, respectively. Employing this ansatz in Eq.(l) we 



arrive at the corresponding Bogoliubov-de Gennes equa- 
tions yielding the excitation spectrum and the Bogoli- 
ubov coefficients Uj and Vj . Interestingly, even in absence 
of hopping, dipolar inter-site interactions result in a col- 
lective character of the excitations that are shared by all 
sites. In consequence, the excitation spectrum acquires 
a band-like character [4] as depicted in Fig. 2. 

Modes with imaginary frequency are associated with 
dynamical instability. For non-dipolar gases, inter-site 
interactions are negligible and hence all transverse modes 
remain degenerated. As a result, modulational instabil- 
ity develops independently in each site and no correlated 
density pattern occurs during the post-instability dynam- 
ics. The situation dramatically changes for sufficiently 
large dipole moment, as the inter-site interactions lift 
the degeneracy between the transverse modes. In partic- 
ular, the most unstable mode becomes significantly more 
unstable than other modes, as shown in Fig. 2, govern- 
ing the BEC dynamics within the linear regime. Cru- 
cially, this most unstable mode is not only characterized 
by a ^-momentum q c (associated with the minimum of 
uj 2 in Fig. 2) setting the modulational instability in each 
wire, but also by a transverse dependence along the y 
direction locking the density pattern between sites. As 
a result, during the first stages of the post-instability 
dynamics a correlated modulational instability develops. 
Interestingly, our numerical simulations predict that this 
phenomenon may be observed in existing Chromium ex- 
periments [16] or even more pronouncedly with recently 
condensed Dysprosium atoms [9]. 

Fig. 3 (top) depicts the case of a 52 Cr BEC desta- 
bilized by an abrupt change of a > into a sufficient 
a < by means of a Feshbach resonance. The numer- 
ical solution of Eq. (1) shows that despite the absence 
of inter-site hopping a correlated density pattern devel- 
ops. As presented in Fig. 3 (top) this instability pat- 
tern survives well into the non-linear regime where the 
density modulation cannot be considered any more as a 
perturbation of the original homogeneous on-site BECs. 
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FIG. 2. (Color online) Bogoliubov spectrum for a 52 Cr BEC 
(/x = 6/is, where \ib is the Bohr magneton) with a density 
cm -3 and a = — 8.5ao (ao is the Bohr radius), occupying 
= 10 sites of a lattice with the inter-site spacing A = 
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512 nm and a lattice depth of 13. 3Er (recoil energy), which 
results in the uj± = 2ty • 26.7 kHz, and l± = 85.3 nm. Here, 
q c = 0.07/l ± . 
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FIG. 3. (Color online) (top) BEC wave function's density 
distribution after 200 ms of time evolution for the same pa- 
rameters as in Fig. 2. For plotting purposes the y- width of the 
tubes has been magnified, (bottom) Dynamics of the Fourier 
transform of the associated column density E(z, £). The dom- 
inating q = peak has been removed for clarity and the re- 
maining distribution has been normalized to the maximum. 
The arrows indicate the harmonics of q c . 



In typical experiments the density alignements may be 
more easily monitored investigating the column density 
^( z ) = ^2m nm ( z )' Contrary to the uncorrected case, 
for which E(z) would show no clear structure, the corre- 
lated instability results in periodically modulated E(z). 
Fig. 3 (bottom) shows the dynamics of the Fourier trans- 
form of E(z,£) which is clearly characterized by the ap- 
pearance of harmonics of q c (compare Fig. 2 and Fig. 3 
(bottom)). 



IV. FILAMENTATION 

The density modulation depicted in Fig. 3 (top) evolves 
into a correlated pattern of solitons. However, the soli- 
tons are created in an excited state, with both internal 
breathing excitation and center-of-mass motion. As a 
result, for insufficient dipolar interactions the correlated 
density modulation is destroyed during the subsequent 
non-linear time evolution. Consequently, the positions of 
solitons at different sites become uncorrelated, not differ- 
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FIG. 4. (Color online) Filamentation of solitons. Here, a 
snapshot of time evolution of the BEC density distribution 
after 500 ms for, in particular, iV m =20 lattice sites, ji= IS/ib 
anda = — 41.7ao- The remaining parameters are chosen as in 
the case of Fig. 2. 



ing qualitatively from the case of non-polar gases. Strong 
inter-site interactions crucially change this picture as the 
correlated solitons in neighboring sites experience an at- 
tractive inter-site potential. Approximating the solitons 
by Gaussians of width <$, such that l ± <C 5 : A, the binding 
energy for two solitons acquires the form 



E b = (-2g d /A 3 )G(6/A), 



(3) 



which differs from the binding energy between point-like 
solitons (— 2#d /A 3 ) by the regularization function 
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with K n the modified Bessel function of second kind. 
As a result of this inter-site soliton attraction, and al- 
though the initial periodicity of the modulation (as that 
of Fig. 3 (top)) is generally lost, self- assembled soliton fil- 
aments form spontaneously (Fig. 4) when the center-of- 
mass kinetic energy of the solitons acquired in the post- 
instability dynamics cannot overcome the binding energy 
given by Eq. (3). 

In order to analyze the dynamical filamentation quan- 
titatively we introduce at this point the time-dependent 
dimer correlation function for sites m and m' 
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and we define the normalized average dimer correlation 
G n (z,t) = G(z,t)/JdzG(z,t), with 
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A proper figure of merit describing the filamentation is 
provided by 



X (t) = G n (0,t)/G n (t), 



(7) 
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FIG. 5. (Color online) (top) Function %(r) for a typical case 
within the fllamentation regime (g = —0.019, ga — 0.0034, 
A = 6 Zjl). In particular, for the parameters that we employed 
in Fig. 2 the time (t = t/cj±) that we here consider equals 
t — 700 ms. (inset) Time-averaged values of \ f° r long times, 
for different values of ga and constant g = —0.019. (bottom) 
Phase diagram of the possible regimes for g < 0, ga > 0. 



where G n (t) = fdzG n (z,t) 2 is the mean value of 
G n (z,t). Such defined function x(£) characterizes the 
tendency of the solitons at different sites to align into a 
filament. 

In the following analysis we consider a simplified case 
of three lattice sites. Fig. 5 (top) shows xif) f° r a typ- 
ical case within the fllamentation regime (see discussion 
below). The sharp initial peak indicates the formation 
of the correlated density pattern shared among all sites 
at the initial stage of the time evolution, as discussed 
in section III. Similarly to Fig. 3, also here the pattern 
is quickly destroyed as the system enters the non-linear 
regime. However, provided sufficently strong dipolar in- 
teractions, the inter-site soliton binding E^ supports the 
formation of soliton filaments and in consequence %(t) 
grows at larger times. Note that x(t) eventually satu- 
rates remaining constant for times typically much longer 
than the usual experimental timescales. 

In contrast, no fllamentation occurs if the dipolar cou- 



pling is insufficient. In this case, at long times xif) aver- 
ages to x — 1 indicating the absence of inter-site soliton- 
soliton correlation. Hence, driving the g^ parameter from 
small to large values results in a transition from a non- 
filamented into a filamented configuration (see the inset 
of Fig. 5 (top)). Ultimately, for a sufficiently large ga the 
repulsive on-site interactions compensate the attractive 
short-range interactions and the system remains stable. 

As a result, we distinguish three distinct regimes of 
dynamics in a stack of ID dipolar gases: (i) unstable un- 
correlated (soliton liquid), (ii) unstable filamented, and 
(hi) stable. As shown in Fig. 5 (bottom), for a fixed value 
of A/l ± these regimes are determined by the ratio ga/\g\- 
For the considered case of three sites and A/l ± = 6 the 
stability boundary line is given by gd/\g\ =0.70, whereas 
the boundary line between the filamented and unstable 
non- filament ated regimes occurs at ga/\g\ = 0.09. For 
the case of 52 Cr (/i = 6{1b) the fllamentation occurs for 
5.2 < \a\/a < 40.2, whereas for 164 Dy (/i = 10/is) it 
occurs for 47.3 < \a\/ao < 367.0. 

We note that for a larger number of sites the system 
is more unstable due to the inter-site attractive interac- 
tions [4]. Also, the boundary between filamented and un- 
stable non-filamented regimes is shifted towards larger g^ 
values due to the enhanced role of the string-like modes 
of the filaments. Hence, even though the qualitative re- 
sults will not be affected, increasing the number of sites 
will in general reduce the fllamentation regime. 



V. CHECKERBOARD SOLITON CRYSTAL 

Interestingly, the sign of g& may be inverted by means 
of transverse magnetic fields [17] or microwave dressing 
in the case of polar molecules [18]. Note that, although 
we consider this case for its theoretical simplicity, quali- 
tatively the same results may be obtained orienting the 
dipoles along the tubes. In both of these cases the emerg- 
ing instability is characterized by the most unstable Bo- 
goliubov mode presenting a staggered ^/-dependence that 
results in an anti-correlated density pattern with maxima 
in a given site aligned with minima in the neighboring 
ones. Strikingly, for a sufficiently strong dipole moment, 
this anti-correlated structure formed at the initial stage 
of the post-instability dynamics seeds the formation of a 
permanent checkerboard soliton crystal in the non-linear 
regime, as shown in Fig. 6. 

Remarkably, while purely repulsive interactions sus- 
tain 2D Wigner-like crystals proposed for polar molecules 
[13, 14], the crystal of solitons is self-maintained by a 
subtle interplay of dipolar inter-tube repulsion and intra- 
tube attraction. Due to the anti-correlated character of 
the density modulation, solitons in neighboring sites pro- 
vide an effective potential barrier that prevents mutually 
attracting solitons in the same tube to come together, 
hence keeping the crystal stable. 

In order to characterize quantitatively the dynam- 
ical formation of a soliton crystal, we employ the 
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FIG. 6. (Color online) Spontaneous crystallization of solitons 
in the case of negative ga- Here, a = 306ao, /h = 36/hb and the 
remaining parameters such as in Fig. 2. 



time windows, as depicted in Fig. 7 (top). For all g^ 
values within the unstable regime the NN anticorrela- 
tion function Xnn^) < 1 remains constant at all times. 
In contrast, depending on the value of the g& parame- 
ter, Xnnn(^) function shows two distinctive types of time 
dependence. Namely, while in the window of the crys- 
tallization regime Xnnn^) saturates at a value indicating 
NNN anticorrelation and so the emergence of a stable 
soliton crystal. Contrastingly, for large dipolar interac- 
tions the initially anticorrelated Xnnn, which originates 
in the linear regime, decreases in time indicating destruc- 
tion of the checkerboard pattern. 

Hence, for negative g& values we identify three dis- 
tinct regimes depicted in Fig. 7 (bottom): (i) a stable 
regime for small dipole values, (ii) an unstable regime 
intrinsically characterized by the dynamical formation of 
a checkerboard soliton crystal, and (iii) a strong dipolar 



notation introduced in section IV, defining the nor- 
malized averaged nearest-neighbor (NN) and next-to- 
nearest-neighbor (NNN) dimer correlations G™(z,t) = 
G a (z,t)/JdzG a (z,t), with a = NN.NNN and 
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and we introduce functions 



Xa{t) = G2(0,t)/G2{t) 



(8) 
(9) 

(10) 

where G£(i) = fdzG™(z,t) 2 is the mean value of 
G£(z,t). The checkerboard soliton arrangement is char- 
acterized by the NN anti-correlation (xmn^) < 1) and 
the NNN correlation (xm\w(£) > 1)- In the following we 
consider a particular case of four lattice sites. A generic 
example of Xnn (t) and Xnnn (t) time evolution within the 
crystalline regime (see discussion below) is depicted in 
the inset of Fig. 7 (top). 

As in the case of filament at ion, the emergence of the 
soliton crystal is limited to a window of \gd\/g values. 
While for a weak dipolar coupling the system remains 
stable, a sufficiently large dipole moment value renders 
the attractive intra-tube interactions dominant and, in 
consequence, we observe the formation of the staggered 
soliton pattern. Note that this configuration, originating 
in the ant i- correlated modulational instability emerging 
within the linear regime, is indeed a highly metastable 
state, as it maximizes NNN dipolar interactions. Cru- 
cially, however, our numerical simulations show that such 
soliton crystal state characterized by Xnn^) < 1 coincid- 
ing with Xnnn^) > 1 remains stable well beyond typical 
experimental timescales, being hence effectively perma- 
nent. Beyond a critical value of the dipolar coupling 
the NNN repulsion destroys the NNN anticorrelation and 
hence the crystal. 

The instability properties of the soliton crystal may 
be studied by considering the average Xa for different 
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FIG. 7. (Color online) (top) The inset shows a typical exam- 
ple of time evolution of Xnn (t) and Xnnn (t) within the soliton 
crystal regime (g = 0.069, gd = -0.032 and A = 6/ x ). We 
average X«(t) within three different time intervals An, 2, 3, 
and we depict in the figure the corresponding averaged x& 
(#, ■, a) for different values of ga and constant g = 0.069. 
Note that for \ga\/g > 0.47, Xnnn decreases in time, indicat- 
ing destruction of the checkerboard crystal. In particular, for 
the parameters we employed in Fig. 2, the time (t — t/uj±) 
that we here consider equals t = 1200 ms. (bottom) Phase 
diagram of the possible regimes for g > 0, gd < 0. 



interactions regime in which only nearest neighbor anti- 
correlation is preserved while the next-to-nearest neigh- 
bor correlation is lost (soliton liquid). In analogy to 
the filament at ion phenomenon, for a fixed A/l ± value 
the regimes boundaries depend solely on the \gd\/g ra- 
tio. For A/Zjl = 6, the crystalliztion regime occurs for 
0.40 < \g d \/g < 0.47, which for 52 Cr ( 164 Dy) requires 
7.7 < a/a < 9.1 (70.0 < a/a < 82.9). 



VI. SUMMARY 

In conclusion, the dipolar inter-site interactions in a 
destabilized dipolar BEC confined in a stack of quasi- 
1D tubes induce an interesting dynamics characterized 
by the development of a correlated modulational insta- 
bility in the non-overlapping sites. For a sufficiently 
large dipole moment this density modulation seeds the 
spontaneous self-assembly of soliton filaments or a soli- 
ton checkerboard crystal, depending on the sign of the 



dipolar interactions. Contrary to filaments and crystals 
of individual molecules, filaments and crystals of solitons 
self-assemble spontaneously merely by simple destabiliza- 
tion of the condensate. Moreover, we expect that due to 
the many-body character of the constituent solitons the 
dipole moment necessary for observing these structures 
may be significantly reduced and that they may be at- 
tainable with partially polarized polar molecules [19] or 
highly magnetic atoms, paving a promising route towards 
the first realization of 2D patterns of solitons in ultra-cold 
gases and, to the best of our knowledge, in nonlinear op- 
tics as well. 
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